import math  
import numpy as np
  
def haversine(lat1, lon1, lat2, lon2):  
    # 将经纬度从度转换为弧度  
    lat1, lon1, lat2, lon2 = map(math.radians, [lat1, lon1, lat2, lon2])  
  
    # Haversine公式  
    dlon = lon2 - lon1   
    dlat = lat2 - lat1   
    a = math.sin(dlat/2)**2 + math.cos(lat1) * math.cos(lat2) * math.sin(dlon/2)**2  
    c = 2 * math.asin(math.sqrt(a))   
    r = 6371 # 地球平均半径，单位：公里  
  
    # 计算结果乘以地球半径  
    return c * r  

# 创建纬度数组  
lat = np.array([  
    40.06694159,  
    31.15012008,  
    30.23345533,  
    34.43349451,  
    26.66676257,  
    24.53353226,  
    25.21674442,  
    18.30014362,  
    22.30025672  
])  
  
# 创建经度数组  
lon = np.array([  
    116.6001157,  
    121.8002472,  
    120.4334292,  
    108.7667046,  
    100.2500782,  
    118.1334968,  
    110.0501363,  
    109.416721,  
    113.9333702  
])  

S = haversine(lat[4],lon[4],lat[3],lon[3])+haversine(lat[3],lon[3],lat[0],lon[0])+haversine(lat[0],lon[0],lat[1],lon[1])+haversine(lat[1],lon[1],lat[2],lon[2])+haversine(lat[2],lon[2],lat[5],lon[5])+haversine(lat[5],lon[5],lat[8],lon[8])+haversine(lat[8],lon[8],lat[7],lon[7])+ haversine(lat[7],lon[7],lat[6],lon[6])+haversine(lat[6],lon[6],lat[4],lon[4])
print(S)